NASA/TM— 1999-209633 



Flutter and Forced Response Analyses of 
Cascades Using a Two-Dimensional 
Linearized Euler Solver 


T.S.R. Reddy and R. Srivastava 

The University of Toledo, Toledo, Ohio 

O. Mehmed 

Glenn Research Center, Cleveland, Ohio 


November 1999 



The NASA STI Program Office ... in Profile 


Since its founding, NASA has been dedicated to 
the advancement of aeronautics and space 
science. The NASA Scientific and Technical 
Information (STI) Program Office plays a key part 
in helping NASA maintain this important role. 

The NASA STI Program Office is operated by 
Langley Research Center, the Lead Center for 
NASA's scientific and technical information. The 
NASA STI Program Office provides access to the 
NASA STI Database, the largest collection of 
aeronautical and space science STI in the world. 
The Program Office is also NASA's institutional 
mechanism for disseminating the results of its 
research and development activities. These results 
are published by NASA in the NASA STI Report 
Series, which includes the following report types: 

• TECHNICAL PUBLICATION. Reports of 
completed research or a major significant 
phase of research that present the results of 
NASA programs and include extensive data 
or theoretical analysis. Includes compilations 
of significant scientific and technical data and 
information deemed to be of continuing 
reference value. NASA's counterpart of peer- 
reviewed formal professional papers but 

has less stringent limitations on manuscript 
length and extent of graphic presentations. 

• TECHNICAL MEMORANDUM. Scientific 
and technical findings that are preliminary or 
of specialized interest, e.g., quick release 
reports, working papers, and bibliographies 
that contain minimal annotation. Does not 
contain extensive analysis. 

• CONTRACTOR REPORT. Scientific and 
technical findings by NASA-sponsored 
contractors and grantees. 


• CONFERENCE PUBLICATION. Collected 
papers from scientific and technical 
conferences, symposia, seminars, or other 
meetings sponsored or cosponsored by 
NASA. 

• SPECIAL PUBLICATION. Scientific, 
technical, or historical information from 
NASA programs, projects, and missions, 
often concerned with subjects having 
substantial public interest. 

• TECHNICAL TRANSLATION. English- 
language translations of foreign scientific 
and technical material pertinent to NASA's 
mission. 

Specialized services that complement the STI 
Program Office's diverse offerings include 
creating custom thesauri, building customized 
data bases, organizing and publishing research 
results . . . even providing videos. 

For more information about the NASA STI 
Program Office, see the following: 

• Access the NASA STI Program Home Page 
at http:/ / www.sti.nasa.gov 

• E-mail your question via the Internet to 
help@sti.nasa.gov 

• Fax your question to the NASA Access 
Help Desk at (301) 621-0134 

• Telephone the NASA Access Help Desk at 
(301) 621-0390 

• Write to: 

NASA Access Help Desk 

NASA Center for AeroSpace Information 

7121 Standard Drive 

Hanover, MD 21076 


NASA/TM— 1999-209633 



Flutter and Forced Response Analyses of 
Cascades Using a Two-Dimensional 
Linearized Euler Solver 


T.S.R. Reddy and R. Srivastava 
The University of Toledo, Toledo, Ohio 

O. Mehmed 

Glenn Research Center, Cleveland, Ohio 


National Aeronautics and 
Space Administration 


Glenn Research Center 


November 1999 



Acknowledgments 


The authors acknowledge the many helpful suggestions from Dr. M. Montogomery and Dr. V. Capece 
during the course of using LINFLX2D. This work is supported by a grant from the 
NASA Glenn Research Center Structural Mechanics and Dynamics Branch under funding from the 
Advanced Subsonic Technology Project and the Turbomachinery and Combustion Technology Project. 
John Rohde, George Stefko, and Kestutis Civinskas are the program managers. 


Trade names or manufacturers' names are used in this report for 
identification only. This usage does not constitute an official 
endorsement, either expressed or implied, by the National 
Aeronautics and Space Administration. 


Available from 


NASA Center for Aerospace Information 
7121 Standard Drive 
Hanover, MD 21076 
Price Code: A03 


National Technical Information Service 
5285 Port Royal Road 
Springfield, VA 22100 
Price Code: A03 




FLUTTER AND FORCED RESPONSE ANALYSES OF CASCADES USING A 
TWO DIMENSIONAL LINEARIZED EULER SOLVER 


T. S. R. Reddy and R. Srivastava 

The University of Toledo 
Toledo, Ohio 43606 

O. Mehmed 

National Aeronautics and Space Administration 
Glenn Research Center 
Cleveland, Ohio 44135 


ABSTRACT 

Flutter and forced response analyses for a cascade of blades in subsonic and transonic flow is 
presented. The structural model for each blade is a typical section with bending and torsion 
degrees of freedom. The unsteady aerodynamic forces due to bending and torsion motions, and 
due to a vortical gust disturbance are obtained by solving unsteady linearized Euler equations. 
The unsteady linearized equations are obtained by linearizing the unsteady non-linear equations 
about the steady flow. The predicted unsteady aerodynamic forces include the effect of steady 
aerodynamic loading due to airfoil shape, thickness and angle of attack. The aeroelastic 
equations are solved in the frequency domain by coupling the unsteady aerodynamic forces to the 
aeroelastic solver MISER. The present unsteady aerodynamic solver showed good correlation 
with published results for both flutter and forced response predictions. Further improvements are 
required to use the unsteady aerodynamic solver in a design cycle. 

INTRODUCTION 

The aeroelastic research program at NASA Glenn Research Center is focused on flutter 
(unstalled, stalled, and whirl), and forced response analyses of propulsion components. An 
overview of this research was presented in Ref. 1. The review showed that a range of 
aerodynamic and structural models have been used to obtain the aeroelastic equations. Both, 
frequency and time domain methods have been used to obtain unsteady aerodynamic forces and 
to solve the aeroelastic equations. It was noted that time domain methods require large 
computational time compared to frequency domain methods, and should only be used when non- 
linearities are expected, or for the final design. 

Two approaches have been used in obtaining the unsteady aerodynamic forces using frequency 
domain methods. In the first approach, Refs. 2-6, the unsteady aerodynamic equations are 
linearized about a uniform steady flow, thereby neglecting the effects of airfoil shape, thickness, 
and angle of attack. The unsteady aerodynamic models developed in Refs. 2-6 have been used in 
Refs. 7-8 to study the flutter and forced response analysis of a compressor cascade using a typical 
section structural model. However, methods developed by this approach are restricted to shock- 
free flows through lightly-loaded blade rows. In the second approach, Ref. 9, the unsteady flow 


* Resident Research Associate, NASA Glenn Research Center, Cleveland, Ohio 
t Research Engineer 


N AS A/TM— 1999-209633 


1 



is regarded as a small amplitude perturbation about a non-uniform steady flow. The unsteady 
non-linear aerodynamic equations are linearized about the non-uniform steady flow, resulting in 
variable coefficient linear unsteady aerodynamic equations, which include the effects of steady 
aerodynamic loading due to airfoil shape, thickness and angle of attack. 

Following the second approach. Refs. 10-11 developed a nonlinear steady and linear unsteady 
aerodynamic model based on the potential equation. This unsteady aerodynamic model was used 
to study the effect of steady aerodynamic loading on flutter stability using a typical section 
structural model in Ref. 12. However, the formulation based on the potential equation requires 
corrections for entropy and flow rotation. The Euler equations can be used to correctly model 
rotational and entropy effects associated with strong shocks. Unsteady linearized Euler 
aerodynamic models that include the effect of steady aerodynamic loading were developed in 
Refs. 13-15. 

Recently, a two dimensional linearized Euler code named LINFLX2D was developed (Ref. 16) 
under a NASA contract. This code is based on the non-linear Euler solver, NPHASE (Ref. 17). 
The objectives of the present study are (1) to couple the LINFLX2D code with the aeroelastic 
code MISER (Refs. 7,8 , (2) validate the aeroelastic predictions, and (3) evaluate the LINFLX2D 
code for use in the aeroelastic design cycle. The MISER code at present calculates the flutter and 
forced response characteristics using the unsteady aerodynamic models based on classical linear 
theories (theories based on flat plate geometry formulation). The MISER code is selected 
because the typical section structural model with pitching and plunging degrees of freedom is the 
basis for aeroelastic formulation. This is in line with the unsteady aerodynamic model of 
LINFLX2D. Over the years, many versions of MISER have been in circulation. To show that 
the version of the MISER code used in the present calculations gives correct results, selected 
cases from Ref. 7 are recalculated and the results presented in Appendix A. 

In the current report flutter and forced response analyses of a cascade of blades is presented using 
the frequency domain approach for a two-dimensional cascade model. Two degrees of freedom, 
plunging and pitching, are considered in the analysis. Vortical gust disturbance is considered tor 
forced response calculations. Brief descriptions of the formulation and method of analysis are 
given in the next section, followed by results and discussion. 

FORMULATION 

The aerodynamic model and the aeroelastic formulation are described in this section. 
Aerodynamic model 

Non-linear Steady Euler Solver, NPHASE 

The steady aerodynamic model is based on the unsteady, two-dimensional Euler equations. The 
equations in conservative differential form are solved in a time-dependent body-fitted curvilinear 
reference frame. This transformation process and the ensuing numerical method are presented in 
detail in Ref. 17. The equations are discretized and solved using a finite volume method using a 
TVD scheme. The steady solutions presented herein are obtained using the implicit scheme 
developed in Ref. 17, which is third order accurate in space and second order accurate in time. 
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Linear Unsteady Euler Solver, LINFLX2D 


To obtain the linearized unsteady Euler equations. Ref. 16, the dependent variables in the 
unsteady non-linear Euler equations are expanded in an asymptotic series of the form 

U = U(x) + u(x(x,t),t) + higher order terms (1) 

where, the term U(x) is of order one and the second term is of the order e. Substituting the 
expansion of Eq. 1 in the nonlinear unsteady Euler equations, and equating terms of like power in 
e, and neglecting terms of second and higher order in e, nonlinear steady equations and linear 
variable coefficient unsteady equations are obtained. The unsteady linear equations are further 
simplified by assuming the unsteady excitations and responses are harmonic in time as 

u(x(x,t),t) = Re[u(x)exp(ia)t) (2) 

For harmonic blade motions with constant phase angle between adjacent blades (interblade phase 
angle), the values of interblade phase angle (a) that can occur are given as (Ref. 18). 

a r = 2nr! N r = 0,1,2 AM (3) 

where N is the number of blades in the cascade. In a time domain approach, the number of blade 
passages required for the solution depends on the interblade phase angle, and small phase angles 
may require a large number of blade passages to calculate the unsteady aerodynamic forces. 
However, with the linear approach, the periodic conditions are applied on a single extended blade 
passage region i.e. a region of angular pitch, 

0 = 2k/N (4) 

In solving the linear unsteady equations, the independent variables are regarded as pseudo time 
dependent. This allows solutions to be determined using conventional time -marching 
algorithms to converge the steady and the complex amplitudes of the unsteady conservation 
variables to their steady state values. For more details, see reference 16. 

Aeroelastic model 

For completeness, the aeroelastic formulation in the frequency domain for a typical section 
structural model (Refs. 7,8) is presented in this section. The analysis is presented for a generally 
mistuned cascade in which each blade may have different structural properties. The analysis for 
the special case of a tuned cascade in which all blades are identical is presented in the subsequent 
section. The approach followed assumes that the structure is vibrating in an aeroelastic mode 
(interblade phase angle mode) with a motion that is a harmonic function of time. The frequency 
of oscillation is permitted to take on complex values thus allowing decaying-, growing- or 
constant-amplitude oscillations. The aerodynamic forces corresponding to constant-amplitude 
harmonic oscillations are inserted into the equations of motion to formulate a complex 
eigenvalue problem. The eigenvalues are generally complex quantities, and therefore a complex 
frequency is obtained. The real part of the complex frequency represents the damping ratio and 
thus its sign determines whether the motion is decaying or growing; the imaginary part 
represents the damped frequency of oscillation. 


NASA/TM — 1 999-209633 


3 



Aeroelastic analysis for mistuned cascade 


The equations of motion for the typical section (see Fig. 1) with structural damping can be 
written in matrix form for the s tb blade as: 

[ M s ]{<? , } + [C, \{q s } = } = {/<?? } = {fas } &) 


where [fcfi } denotes motion-dependent aerodynamic loads and {f a } denotes motion- 
independent aerodynamic loads. 


For the two degrees of freedom considered here, equation 5 can be rewritten as 


1 *«s 
*«s rls 



+ 


^tyisChs 0 j b s lb \ 
0 2rasGhsCas. \ Os I 


alg 0 \h s lb \ _ f f hs im s b \ 

0 ^as<U^s_ \ I \faslm s b^ j (6) 


where h is the plunging (bending) displacement normal to the chord, a is the pitching (torsion) 
displacement, xa the distance between the elastic axis and center of mass in semi-chord units; r a 

is the radius of gyration about the elastic axis in semi-chord units; £/j and £« are the damping 
ratios; b is the airfoil semi-chord; fit)/, is the uncoupled natural frequency for bending; tu« is the 
uncoupled natural frequency for torsion; fh and fa are the aerodynamic loads. 

It is assumed that the motion of the blades is harmonic in time with a frequency (0 and is 

% 

given by 



Note that the motion has been represented as the sum of contributions from each interblade phase 
angle mode in which each blade has an amplitude h a ^b, 0Ca r and the phase angle between 
adjacent blades is given by Eq. 3. 

The corresponding aerodynamic forces can be written in terms of the complex-valued unsteady 
aerodynamic coefficients // ih , la} i , lha . lota » Iwh* ar| d Iwtt ( see Ref- 7 for the definition of 
these coefficients) as 

fhs t m s b 1 _ py_ 'ST l \jhhrKr ! b + ^ar^ar]] , Ua^a.s + ® | [C*r]l ita^a.s ^ 

LJ m f) Us hs \[Lr h ar /b + l aar a ar]\ [[far]) 

where fis = ms I Kp„b is the mass ratio of the blade, and p <*> is the density of air. 
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Using Eq. 7 and Eq. 8, Eq. 6 can be written as 


~[M s ][ hoslb \e Mt + A [K s ] ! A °* /6 ] c «* 
I « 0S / L «o.s / 

= I [/^ | + X 

r=0 ' / r =0 

where 


[M s ] — /i s 


1 x a 

[«,]= 

(°hs^ c 1) (1+2 iQts) 

0 

2 

x a r a_ 


0 

(<%/<$) (1+2 it^g) 


[Ar] = 


^/zAr hi ar 
^ahr l(xar. 


^ADA=^ l Wkr ^ 

\ r , \l ‘ 


v w ar 


2 

A = (tub /(o) 


where (o 0 is the reference frequency and the damping terms are approximated as, 


2 i (tiCOhs Chs ~ 2i G)fix i^/is 
2 if asQXOfxs^as ~ 2iv £ s 0)(xs^ t 


as 


(10) 


To proceed further, the equations for all the N blades on the disk must be considered. For the 
assumed harmonic blade motion, the displacements ( X } of all blades can be written as a sum of 
contributions from all interblade phase angles as 


{X)e im = [E}{Y}e“* . (11) 


where { Y } consists of the displacement amplitudes corresponding to the rth interblade phase 
angle and [E] is the transfer matrix, the elements of which are given by 

E(s,r) = e 2nsri /N 


where s is the blade counter and r is the interblade phase angle counter. 
Using this relation, we obtain: 


-[ £ fMPl!i'! + a [Bp [*][£] (r i = [A]jy | + [ad\ 

where 


( 12 ) 


[M] is the mass matrix for all blades with diagonal elements consisting of mass matrix of each 
blade. Similarly [K] is the stiffness matrix and [A] is the aerodynamic matrix for all blades, 
{AD} is the aerodynamic forces due to wake on all the blades. 
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Finally, after rearranging, the equations can be written as: 

[[P]-A[Q]](r}“-{AD) (13) 

where 

[P] = [[£]-' [M][£] + [A]] 

[Q] = [£]■' M[B] 

For a stability calculation (flutter), the motion-independent forces {AD} are set to zero and the 
eigenvalue problem is obtained in the standard form: 

[[/>]- a[ 0 ]Rk}={o} < 14 > 

The solution of the above eigenvalue problem Eq. 14 results in 2 N complex eigenvalues of 
the form 


iSH. - — L— = /j + i v 

"o VT (15) 

The real part of the eigenvalue ( Ji ) represents the damping ratio, and the imaginary part ( v ) 
represents the damped frequency; flutter occurs if ^ 0 for any of the eigenvalues. 

The aeroelastic response of the blades induced by wakes is calculated from equation (13) as 

jy) = -[[P]-A [<?]]*>!>) 06 ) 


The amplitude of each blade is obtained by substituting equation (16) into equation (11). 


Aeroelastic analysis for a tuned cascade 


For a tuned cascade (or rotor), in which all the blades are identical, the foregoing analysis can be 
simplified considerably. In this case, the aeroelastic modes consist of individual blades vibrating 
with equal amplitudes with a fixed interblade phase angle between adjacent blades. Hence, for 
this problem, the motion of the typical blade is written as 


thjb) \K 


a. 



i(tit i<J J 

e e 


(18) 


Thus the equation for the blade becomes 



+ 


*[*■]{! 'C 


-M 


\K,' b \ 

a... 


gH<tl+G r .s) 


+ 


{AD r }e 


i( ftJ + G r s ) 


(19) 
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Since the blades are identical, the same equation is obtained for each blade. Thus, no additional 
information can be obtained by assembling the equations for all the blades on the disk as was 
done for the general mistuned system. Instead, equation (19) is solved for N different values of 
the interblade phase angle given by equation (3). As before, the equations for the stability 
(flutter) calculation are obtained by setting the motion-independent forces to zero. 

For the stability calculation, the equation can be simplified as 

(20) 

where 


[Pr] = 


where the subscript 's' identifying the blade has been dropped and the reference frequency (0 o has 
been chosen to be equal to the torsional frequency a>cc . 

The solution of the above eigenvalue problem results in two complex eigenvalues of the form 
fi + i v , and flutter occurs if fl > 0 . For the tuned cascade, the stability of each phase angle 
mode is examined separately. Hence, the interblade phase angle is fixed at one of the values 
given by equation (3), and the 2X2 eigenvalue problem is solved. The value of interblade phase 
angle is then changed, and the procedure is repeated for each of the N permissible values. The 
critical phase angle is identified as the one which results in the lowest flutter speed. 

Stability calculation 

The aerodynamic coefficients are calculated before the eigenvalue problem can be set up and 
solved. Since the unsteady aerodynamic coefficients depend on the frequency of oscillation, it is 
necessary to assume a frequency 0) (reduced frequency of blade vibration based on chord, k c = 

<a c / U, where U is the free stream velocity ) in advance to be able to calculate the aerodynamic 
coefficients. In actual calculations, the aerodynamic coefficients are functions of inlet Mach 
number Moo , and interblade phase angle o>- , in addition to cascade geometric parameters. In the 
current study, for a given mlet Mach number, the reduced frequency is varied until the real part 
of one of the eigenvalues fl becomes zero while the real parts of the remaining eigenvalues are 
negative. The assumed flutter-reduced frequency k c f and the calculated flutter frequency Vf 
are both based on oyf . Thus, these two can be combined to eliminate <*yand the flutter speed is 

obtained, namely, Vf = VfC &b / k c f . Since the inlet Mach number is known, this flutter speed 
gives the inlet condition (speed of sound, aoo) at which the cascade will be neutrally stable for the 
given Mach number. This procedure can be repeated to obtain a plot of flutter speed versus 
Mach number. Knowing the operating conditions, it is possible to determine whether flutter will 


I U~* r l'hhr 

^(ft) ; /ftfc) 2 ( l+2i&) 
l cx~^~ lothr 

Hrh ( 1 + 2 iCa) 


fKo^/aJ (\+7iQ 

fir a (l+2i£*) 


NASA/TM— 1999-209633 


7 



occur within the operating region and it so, the Mach number and frequency at flutter. It should 
be noted that the analysis of a mistuned cascade requires the solution of the equations for all 
phase angles at one time for a given reduced frequency. 

RESULTS AND DISCUSSION 

The analysis methods presented earlier are now applied to investigate the behavior of cascades 
oscillating in pitch and plunge, and cascades subjected to unsteady vortical gust disturbance. 
First, unsteady pressure distributions and loading for selected cascade configurations are 
presented and compared with published results to help validate the present linearized unsteady 
Euler solver, LINFLX2D. Then, flutter and forced response predictions with and without 
mistuning are presented. For the calculation of unsteady aerodynamic forces, first NPHASE, the 
non-linear Euler solver is used to generate the steady solution. This steady solution is then used 
as an input to LINFLX2D, to calculate the unsteady aerodynamic forces due to small unsteady 
perturbations in pitch, plunge and vortical gust disturbance. The calculated unsteady 
aerodynamic coefficients are used in MISER to predict the flutter stability and forced response. 
The MISER code used in the present investigation has been verified and sample results are 
presented in Appendix A. 

Code Validation 

Unsteady pressure difference, the difference between the lower and upper surfaces of the airfoil, 
were calculated for a flat plate cascade, and a cascade designated as the tenth standard 
configuration (CIO), Ref. 16. The unsteady pressures difference was normalized by ((air 
density* square of free stream velocity*abs(amplitude of pitching motion + amplitude of 
plunging motion* kc ))• Both cascades have a stagger angle (0) of 45 degrees, and a gap to chord 
ratio (s/c) of 1.0. The tenth standard configuration airfoils are constructed by superposing the 
thickness distribution of a modified NACA 5506 airfoil on a circular-arc camber line. See Ref. 
16 for more details. Two Mach numbers are considered, M=0.7 and M=0.8 for a reduced 
frequency of oscillation based on semichord (k h ; k b = 0.5* k c ) of 0.5. At M=0.7 the flow is 
subsonic, and at M=0.8 the flow is transonic with a shock at 25% chord from the leading edge of 
the airfoil. The steady angle of attack is zero for the flat plate geometry. It is ten degrees for the 
CIO airfoil at M=0.7, and thirteen degrees at M=0.8. Unless otherwise stated, a 141x41 H-grid is 
used for the study. There are 80 points on the airfoil, the inlet boundary is 5 chords from the 
airfoil leading edge, and the exit boundary is 10 chords from the airfoil trailing edge. Some 
calculations were also done with a 155x41 grid for which the inlet and exit boundaries are 
located at one chord from the leading and trailing edges, respectively. There are 55 points on the 
airfoil for this grid. 
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Unsteady Aerodynamic Pressures 


Subsonic Inflow 


Figure 2a shows the unsteady pressure difference distribution obtained for the flat plate cascade 
operating at M = 0.7. The interblade phase angle, a, is 180 degrees, the reduced frequency of 
oscillation (k b ) is 0.5, and the Mach number is 0.7. The blades are oscillating in pitch about the 
midchord. For the inlet Mach number of 0.7, the flow is shock free. Figure 2a shows predictions 
form linear theory (Ref. 3) and from the present linearized Euler code, LINFLX2D. The 
predictions from LINFLX2D correlate well with the linear theory results. 

Figure 2b shows the unsteady pressure difference distribution obtained for the standard cascade 
configuration, CIO. The flow to the cascade is at 10 degrees angle of attack. Again, the 
interblade phase angle is 180 degrees, the reduced frequency of oscillation (k h ) is 0.5, and the 
Mach number of 0.7. The blades are oscillating in pitch about the midchord. For the inlet Mach 
number of 0.7, the flow is shock free. Figure 2b shows predictions from linear theory (Ref. 3), 
nonlinear Euler (Ref. 17) and from the present linearized Euler code, LINFLX2D. The 
predictions from LINFLX2D correlate well with the nonlinear Euler results. As expected, there 
is a difference with linear theory, since the effects of airfoil geometry and angle of attack are not 
included in the linear theory. 

Transonic Inflow 


For the standard configuration, CIO, and for an inlet Mach number of 0.8, the flow is transonic 
with a normal shock occurring in each blade passage. The steady angle of attack is 13 degrees. 
The steady Mach number distribution for this flow obtained from NPHASE is shown in Fig. 3a. 
There is a normal shock at about 28% of the chord. Results from the full potential solver. 
Ref. 19, are also included for comparison. Both results agree well, except that the nonlinear 
steady Euler solver, NPHASE, predicts the shock location slightly downstream of that predicted 
by the full potential solver. A shock fitting procedure was used in Ref. 19, whereas the shock is 
captured in the present solver. 

The unsteady pressure distribution from LINFLX2D is shown in Fig. 3b, along with a 
comparison with linear theory (Ref. 3) and from the unsteady nonlinear Euler solver (Ref. 17). 
The unsteady results are for pitching about mid-chord with kb = 0.5, <7 = 180° and amplitude of 
oscillation,^ = 2°. As expected, linear theory does not show any shock, indicating that 
unsteady analysis based on linear theory will not be accurate for cascade flutter analysis in 
transonic flow. Some discrepancy is seen between the unsteady results obtained from 
LINFLX2D and from the unsteady nonlinear Euler solver. For the transonic flow case the 
following factors may effect the results more than for subsonic case: (1) grids used, (2) accuracy 
lost when time series data is processed to obtain frequency data or (3) the steady solution on 
which LINFLX2D based is not a low loss solution. Even though these issues need separate 
study. Fig. 3b shows a fair correlation between LIFLX2D results and non-linear unsteady Euler. 

Similar results were obtained for plunging motion for both subsonic and transonic inflows. 
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Unsteady Aerodynamic Force coefficients 

It is the unsteady aerodynamic force coefficients, unsteady lift and moment that are obtained by 
integrating unsteady pressure differences, which are used in the flutter and forced response 
prediction. Table 1 shows the unsteady aerodynamic force coefficients obtained using the Smith 
code, Ref. 3, the present LINFLX2D code, and the non-linear Euler solver of Ref. 17. In table 1, 
the first column shows the geometry of the airfoil studied, the fourth, fifth and sixth columns 
show the force components for pitching motion, plunging motion and for vortical gust 
respectively. The symbols in the second column are cl for lift coefficient, and cm for moment 
coefficient, with moment taken about the mid chord (x„ = 0.5). The rows show the geometry and 
the Mach number studied. The analysis method is given in column 3 with the notation defined at 
the bottom of the table. The geometry and unsteady conditions are given on the top of the table. 

A comparison of the force coefficients for flat plate geometry in the first row shows that the 
prediction by LINFLX2D are very close to those given by Smith code. For a flat plate geometry 
in subsonic flow, the Smith code gives the exact answer. It is to be noted that the flat plate is 
unloaded and shock free. The values of the force coefficients in the second row show some 
difference from the Smith code, since in LINFLX2D the effect of airfoil and angle of attack are 
included. The force coefficients in the third row include the effect of shock and its motion in 
addition to the effects of airfoil shape and angle of attack. The force coefficients obtained with 
LINFLX2D code compare well with those obtained from a non-linear unsteady Euler solution, 
indicating that LINFLX2D code can be used to analyze cascades in transonic flow. It is to be 
noted that only limited runs for vortical gust response predictions have been made. It is also to 
be noted that the solution for transonic flow could not be obtained with 155x41 grid considered 
here for which the inlet and exit boundaries are at one chord length. It was suggested that for the 
LINFLX2D code to work properly, it is better if the inlet and exit boundaries are at least 
5 chords from the leading and trailing edges respectively (Ref. 20). 

Flutter Calculations 

As mentioned before, given the unsteady aerodynamic force coefficients, the MISER code 
(Ref. 7) calculates the aeroelastic stability of tuned and mistuned cascades oscillating in pitch 
and plunge. The aeroelastic stability results from MISER are presented as a root locus plot, with 
real part of the eigenvalue indicating damping and the imaginary value indicating frequency for 
all possible interblade phase angles. In the present study the unsteady aerodynamic force 
coefficients obtained from LINFLX2D are coupled with MISER as follows. First a steady 
solution is obtained from NPHASE for the required cascade geometry and flow conditions. With 
the steady solution as input, LINFLX2D is run for given a reduced frequency, interblade phase 
angle, and mode of interest (pitching or plunging). The unsteady aerodynamic coefficients 
obtained from LINFLX2D are then stored in a data base, and later used in MISER to calculate 
the stability characteristics. 

The structural properties are the mass ratio, (X, is 258.50; ratio of bending to torsion 
frequencies , (Oh /co a is 0.357; the radius of gyration about the elastic axis in semi-chord units, r a , 
is 0.5774; with no structural damping. The elastic axis (see Fig. 1) is at midchord, i.e. ah = 0.0. 
The distance between the elastic axis and center gravity, xqc, is zero, indicating that coupling 
between bending and torsion is very weak and the flutter mode is dominated by torsional motion. 
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Therefore, only torsion mode is shown for all the calculations. It should be noted that these 
structural properties were also used in Refs. 7 and 8. 

The following procedure is used for stability prediction. First MISER is run for the selected 
cascade (rotor) using linear theory ( Ref. 3) to identify the flutter point. Then MISER is run 
using the LINFLX2D data base for this flutter condition. The eigenvalues obtained from linear 
theory and LINFLX2D are then compared. 

Flat Plate Cascade 


To validate the coupling of the unsteady aerodynamic coefficients obtained from LINFLX2D 
with MISER, a nine bladed cascade with a flat plate geometry was selected. The flutter results 
obtained are compared with those obtained from linear theory. The gap-to-chord ratio is unity 
and the stagger angle is 45 degrees. For the nine bladed cascade, for a given reduced frequency, 
the unsteady force coefficients have to be calculated for the possible interblade phase angles of 0 
40, 80, 1 20, 1 60, 200, 240, 280, and 320 degrees. 

Using linear theory, for M = 0.7 with elastic axis at midchord, for the nine bladed cascade 
mentioned above, an instability was found at k h = 0.18 and at a phase angle of 80 degrees. 
Therefore, LINFLX2D was run for this value of reduced frequency, and for all allowable phase 
angles, and for both plunging and pitching oscillations, and a database was prepared. Then 
MISER was run to predict the stability boundary with the unsteady aerodynamic force coefficient 
data base obtained from LINFLX2D. 


The root locus plot showing the real and imaginary parts of the eigenvalues for the flat plate 
geometry obtained for torsion motion is given in Fig. 4. The unsteady force coefficients for eight 
phase angles were obtained with LINFLX2D. For the o’ = 40 degrees, lor which acoustic 
resonance was expected, the LINFLX2D code did not give a solution for the grid, oscillation of 
amplitude, and time step considered. This case will be studied further in future. It was noted 
that for a flat plate geometry, the computational times were higher than those for airfoil geometry 
(presented below). Figure 4 shows the comparison between predictions from linear theory and 
from LINFLX2D. It can be seen from Fig. 4 that LINFLX2D compares very well with linear 
theory results for the flat plate geometry. The plot shows that the cascade is unstable at the 80' 
phase angle mode. 

Mistuned Cascade 


Results are presented now for cascades with mistuning. The type of mistuning considered is the 
one in which the odd and even numbered blades have different torsional frequencies, also known 
as alternate blade mistuning. This has been studied in Refs. 7 and 8 with 1% and 1.5% alternate 
mistuning. For 1% mistuning, the frequency ratio co a / 0 ) (I is 1.005 for all the even blades, and is 
0.995 for all the odd blades. The study in Ref. 7 and 8 considered a 56 bladed cascade, and 
showed that alternate mistuning has a considerable effect on stability. A 12 bladed cascade is 
considered for the present study, to limit the computational time. MISER was run first for a 
tuned cascade with linear theory, for M=0.7. It was found that the cascade is unstable for k b = 
0. 1 and a = 60 . This value of k h = 0. 1 is used for all the calculation presented below. It should 
be noted here that for this cascade acoustic resonance was expected for o = 0" and 30“. However, 
for the CIO airfoil geometry unsteady solutions were obtained with LINFLX2D for these 
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interblade phase angles also, unlike for the flat plate geometry for which solutions could not be 
obtained with LINFLX2D at acoustic resonance (as mentioned in the previous section). 


First, the amount of mistuning required to change the root locus plot is investigated, i.e. whether 
the effects of 1% and 1.5% alternate as reported in Ref. 7 also occurs for this cascade. It should 
be noted here that the results from the present MISER code were checked with those given in 
Ref.7 for the 56 bladed cascade before starting the current calculations, and are given in 
Appendix A for both flutter and forced response calculations. Figure 5 shows the root locus plot 
with 1%, 1.5%, 5%, 10% and 20% alternate mistuning for the flat plate geometry tor M - 0. / an 
k = 0.18. It can be seen that an alternate mistuning of 1-10% did not show any effect on the 
eigenvalues. Only 20% alternate mistuning showed some effect on the eigenvalues. This is in 
contrast to the level of mistuning used in Ref. 7. However, Ref. 7 has considered a 56 bladed 
cascade in incompressible flow for which the cascading effect is much stronger. 


However, the large amount of mistuning required here to have a significant effect on the 
eigenvalues (shown in Fig. 5) is not surprising. Similar values of mistuning have been used in 
Ref 8 That study used an amount of 1 to 20% mistuning. An examination of Fig. 1 1 ot Re . 
showed that the effect of alternate mistuning decreased with increasing Mach number. In the 
context of Ref. 8 and the present results, it is concluded that for the Mach numbers considere 
here a mistuning of 20% or more has to be used to show a significant effect of on the 
eigenvalues. Therefore, the rest of the calculations are done with 20% alternate mistuning. 

Figure 6 shows the root locus plot obtained for the CIO airfoil for both tuned and 20% alternate 
mistuning, for M = 0.7 and k b = 0.1 . Comparison with flat plate geometry is also shown. The 
following are noted from Fig. 6: (1) For the tuned cascade both flat plate geometry cascade and 
CIO airfoil cascade are unstable for a = 60". However, the CIO cascade is less unstable than the 
flat plate cascade. This indicates that for the tuned cascade the airfoil steady loading effects are 
stabilizing, compared to that for a flat plate airfoil; (2) the addition of 20% alternate mistuning 
moved the cascade with CIO geometry from an unstable to a stable condition, but the flat plate 
cascade remained unstable; (3) the root locus plot for the CIO airfoil with mistuning has two 
frequencies identified with a a = 240" mode, whereas, the O = 60" mode is completely missing. 
The reason for this behavior is under investigation. 


The root locus plot for M = 0.8 and k„ = 0.1 is shown in Fig. 7. The following are noted form 
Fig. 7: (l)the tuned flat plate geometry cascade was unstable for o = 90' and the CIO airfoil 
cascade was stable. This indicates that the airfoil steady loading effects are stabilizing for 
the tuned cascade, compared to that for a flat plate airfoil; (2) the addition of alternate 
mistuning made the CIO geometry more stable, but the flat plate cascade remained unstable; 
(3) For the C 10 airfoil it was noted that for tuned case both a = 90° and G = 120 modes have the 
same amount of damping, but with mistuning, a = 90" mode became more stable than the 

a = 120" mode. 


Forced Response Calculations 

In the formulation for response prediction in MISER, it is possible to consider an excitation 
function consisting of all harmonics of rotational speed of the rotor which range up to r - N- 1 . In 
engine aeroelastic terminology, the harmonic number r is known as the “ engine order of the 
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excitation. The number of harmonics and interblade phase angles are related by Eq. 3 as 
Or = 2nnN r = o, ],2, ••• , AT— i . The coefficients l whr and l WOr in Eq. 8 represent the forcing 
functions in the bending and torsion equations, respectively. To understand the nature of the 
response, excitation in only one harmonic at a time is considered. This results in no loss of 
generality since the principle of superposition holds. If the r = R harmonic is considered, then 
the column matrices {ADJ, { AD, } { AD N _, } are zero except {AD r }. Generally, a limited 
number of engine order of excitations are of interest in forced response calculations. Thus, all 
harmonics or all interblade phase angle modes need not be considered in forced response 
analysis, as in flutter analysis. In the following only torsional response amplitudes are shown 
since bending amplitudes are very small in the range of excitation frequency considered. 

Figure 8 shows the tuned aeroelastic response for R = 6 i.e. a = 180 degrees for flat plate 
geometry. Again a twelve blade (N = 12) cascade with the same structural properties used for 
flutter calculations is considered. However, a structural damping ratio of 0.002 is included in the 
analysis. The unsteady aerodynamic coefficients are obtained at kb = 0.5 for which the cascade 
is aeroelastically stable in all the modes. Calculations are done with a 155x41 grid for which the 
inlet and exit boundaries are located at one chord from the leading and trailing edges 
respectively. There are 55 points on the airfoil for this grid. The forcing frequency range 
investigated is limited to a small range around the uncoupled torsional frequency. For the tuned 
cascade, the response will be entirely in the r = R mode, and all the blades will have equal 
amplitudes. The amplitude of response is a function of to / 0) o . Figure 8 shows the torsional 
response obtained using linear (Smith) theory and present LINFLX2D code for vortical gust 
disturbance. The moments are taken about leading edge (ah = -1.0). Two responses are shown 
in Fig. 8 one obtained by including the motion dependent aerodynamic matrix [A] in Eq. 12 and 
the other without the contribution from [A], The terms of [A] provide aerodynamic damping in 
response calculations. It can be seen that calculations from Smith theory and from LINFLX2D 
are identical indicating that the coupling of LINFLX2D coefficients to MISER code is accurate. 

The inclusion of [A] in response calculations requires the calculation of motion dependent 
unsteady aerodynamic coefficients for all interblade phase angles. However, as mentioned 
earlier, in forced response calculations, only a limited number of interblade phase angles (engine 
orders) are of interest and the calculations of motion dependent unsteady aerodynamic 
coefficients is unnecessary. This essentially ignores the effect of aerodynamic damping on 
response. As forced response is usually calculated for an aeroelastically stable system (positive 
aerodynamic damping), this gives a conservative measure on the response calculations. 
Therefore for further calculations the elements of matrix [A] are not included in the calculation. 

Figure 9 shows the torsional response obtained for the standard configuration in subsonic flow 
with and without 20% alternate mistuning. Comparing the response of CIO with that of flat 
plate, it can be seen that the response has increased considerably. This may be due to high steady 
loading on the airfoil. Also shown in Fig. 9 is the response with 20% alternate mistuning. The 
response has split into two resonance peaks, one for odd blades and the other for even blades, 
resonating at the 0.9 and 1.1 frequency ratios, as expected with no difference in amplitude. It 
should be noted that the study of Ref. 7 showed that when the aerodynamic damping terms are 
included in the response calculations, the amplitudes could be lower or higher depending on the 
engine order. 
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The calculations given above for flat plate geometry show that the LINFLX2D data base is 
accurately coupled with MISER for forced response calculations, as was shown for flutter 
calculations in the previous section. 

CONCLUDING REMARKS 

A two-dimensional unsteady linearized aerodynamic Euler solver, LINFLX2D was used to 
calculate unsteady aerodynamic forces on oscillating cascades. The unsteady aerodynamic forces 
were used in the aeroelastic code MISER to calculate flutter and forced response characteristics 
of cascades. Results were presented for cascades in subsonic and transonic flow. Calculations 
with the linearized Euler solver showed good correlation with the published data. The study also 
showed that the unsteady aerodynamic force data base (for flutter and forced response) from 
LINFLX2D is successfully coupled with MISER to calculate flutter and forced response of tuned 
and mistuned cascades. The following were observed during the study. 

(1) The airfoil steady loading stabilized the cascade for the Mach numbers considered. 

(2) A relatively high amount of alternate mistuning is required to have an effect on the stability 
of the cascade for the Mach numbers considered. 

(3) The steady loading due to the airfoil shape and angle of attack increased the blade response 
compared to that for a flat plate. 

(4) An accurate steady solution called “low loss solution “ is required for the unsteady solver 
LINFLX2D to work. This requires a considerable amount of computational time. This may 
be the biggest drawback of using this code in a design cycle. 

(5) The success of getting an unsteady solution from LINFLX2D also depended on the 
grids used. 

(6) For predicting flutter boundary, calculations have to be done for a range of reduced 
frequencies. It is suggested that LINFLX2D data base be prepared for three or four reduced 
frequencies, and then use interpolation for flutter calculations. These reduced frequency 
values can be selected by first running MISER with linear theory. 

(7) The unsteady aerodynamic solution times varied from 20 minutes to 75 minutes on an SGI 
workstation. This is directly related to the number of iterations required for convergence. A 
solution for flat plate geometry took more time than for the CIO airfoil geometry. 

(8) The LINFLX2D code can be used with any aeroelastic code that uses a typical section 
aeroelastic model as its basis. It can be used with the NASA Glenn Research Center codes 
ASTROP2, and FREPS with little modification. 
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APPENDIX A 


Calculated results from the current Version of MISER 

Over the years, many versions of MISER code have been in circulation. The present authors also 
made some modifications to MISER. An important one was updating calls for numerical 
routines from IMSL 9 to IMSL 10 and LAPACK numerical packages. So, calculations are 
shown in this appendix to confirm that the authors’ MISER code gives the same results as the 
MISER code used in Ref. 7 in 1982. 

Flutter and forced response results are reproduced here by running the current version of the 
MISER code for the example cases given in Ref. 7. A 56 blade rotor representative of a 
compressor stage is used. The stagger angle is 54.4 degrees and the gap to chord ratio (s/c) is 
0.534. The structural properties are the mass ratio, p, is 258.50, ratio of bending to torsion 
frequencies a>b /a>a is 0.357, the radius of gyration about the elastic axis in semi-chord units > a 
is 0.5774, with no structural damping. The elastic axis is at midchord i.e. ah = 0.0. The distance 
between the elastic axis and center gravity, xa is zero. Thus coupling between bending and 
torsion is very weak, and the flutter mode is dominated by torsional motion. Therefore, only the 
torsion mode is shown for all the calculations. The flow is considered incompressible. 

Figure A.l shows the root locus plot obtained with and without mistuning for kb = 0.6425, and 
for zero structural damping. These results are identical with Fig. 4 in Ref .7. 

In figures A2 and A3, aeroelastic response for two engine orders, R = 1 1 and 39, are presented. 
These engine orders correspond to interblade phase angles of 70.7 and 250.7 degrees 
respectively. An alternate mistuning of 1%, and a structural damping of 0.002 is used in the 
calculations, with kb = 1.2. These figures are identical with those shown as Fig. 11a and 1 lb 
of Ref. 7. 
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Table 1 

Unsteady Aerodynamic Force Coefficients 


s/c = 1.0, stagger = 45", a = 180°, k„ = 0.5, pitching about midchord (x0=0.5); steady angle of attack 
is zero for flat plate, 10 degrees for M=0.7 and 13 degrees for M=0.8 for C 10 airfoil. 


Geometry / Flow 
condition 

Force 

Coefficient 

code 

Plunging 
Real, imaginary 

Pitching 
Real, Imaginary 

Vortical gust 
Real, Imaginary 


cl 

a 

-3.5675,1.2945 

-4.2440,0.735 1 

1.5768,-2.8183 

Flat plate 


b 

-3.5133,1.2338 

-4.1786,0.6609 

1.5684,-2.7637 

cascade, M= 0.7 







— 

cm 

a 

0.5595,-0.6765 

0.7138,-0.8777 

-0.1496,0.5944 



b 

0.5297,-0.6417 

0.6780,-0.8324 

-0.1409,0.5648 



■ 

-3.5675,1.2946 

-4.2439,0.7351 

1.5768,-2.8183 

CIO cascade 

cl 

lit 

-2.5502,-0.0132 

-2.2566,-0.1537 

4.7945,-4.6932 

M = 0.7 


Ip': 

-2.6370,0.0553 

-2.3040,-0.0487 

** 



■ 

-2.6576,0.0310 

-2.4415,-0.0362 

** 



■ 

0.5595,-0.6765 

0.7138,-0.8777 

-0.1496,0.5944 



H 

0.6497,-0.3330 

0.6001,-0.4540 

-0.5963,0.4450 


cm 

m 

0.6124,-0.3340 

0.5513,-0.4406 

** 



19 

0.6346,-0.3474 

0.5917,-0.5108 

** 



- 

a 

-0.4792,2.6405 

-0.4467,3.1039 

-0.5932,-2,1091 

CIO cascade 


b 

* 

* 

* 

OO 

o 

ll 

£ 

cl 

c 

-2.5945,1.0469 

-1.8953,0.8885 

** 



d 

-2.7240,1.1272 

-2.0532,0.9119 

** 



a 

-0.4633,0.03245 

-0.7879,-0.1658 

0.02522,-0,2289 



b 

* 

* 

* 


cm 

c 

0.3352,-0.6922 

0.2215,-0.6712 

** 



d 

0.3615,-0.6567 

0.2657,-0.6940 

** 

a = linear theory, 1 

ilef. 3; 


b = LINFLX2D (155x41 grid with 55 points on the airfoil, inlet=-1.0, exit=+2.0) 
c = LINFLX2D (141x41 grid with 81 points on the airfoil, inlet=-5.0, exit=+ 10.00) 
d = nonlinear unsteady Euler code, Ref. 17 

* could not get solution with 155x41 grid; ** did not run for gust predictions 
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Figure 1 : Typical section blade model with two degrees-of-freedom. 



distance along the chord, x/c 


Figure 2a: Unsteady pressure difference distribution for a flat plate geometry cascade; pitching 
about mid-chord; gap-to-chord ratio, sic = 1 .0, stagger angle, 6 = 45°, free stream Mach number, 
M = 0.7, steady angle of attack, / = 0°, reduced frequency based on semi-chord, kfj = 0.5, 
interblade phase angle, <7= 180°, amplitude of oscillation, (Xq = 2°. 
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Figure 2b: Unsteady pressure difference distribution for a subsonic cascade, pitching about 
midchord; , s/c® 1.0, 0=45°, AT =0.7, k b = 0.5, <7= 180°, 0 ^ = 2°, / = 0° for flat plate and 

i= 10° for CIO airfoil. 



Figure 3a: Steady Mach number distribution for a transonic cascade; standard configuration 
No. 10,5/c= 1.0, 0=45°,M = O.8, i = 13°. 
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Figure 3b: Unsteady pressure difference distribution for a transonic cascade; standard 
configuration No. 10, pitching about mid-chord; parameters as in Fig. 3a, kb = 0.5, cr= 180°, 
= 2 °. 



real part of eigenvalue 


Figure 4: Root locus plot for 9 bladed cascade, torsional motion, flat plate geometry, kfy = 0.18, 
a h = 0.0, )i = 258.5, r a = 0.5774, structural damping=0.0, (Of , /(O a = 0.357, M = 0.7, / = 0° 
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. 0.2 - 0.15 - 0.1 - 0.05 0 0.05 0.1 

real part of eigenvalue 

Figure 5: Effect of alternate mistuning on eigenvalues in subsonic flow; flat plate geometry, 
M = 0.7, i = 0°, kb - 0.1, N=12, pitching about mid-chord; .1%,2%,5%,10%, 20%; structural 

parameters same as in Fig. 4. 



- 0.2 - 0.15 - 0.1 - 0.05 0 0.05 0.1 

real part of eigenvalue 

Figure 6: Effect of alternate mistuning on eigenvalues in subsonic flow; standard configuration 
No. 10, M = 0.7, i = 10°, kb~ 0.1, N = 12, pitching about mid-chord; structural parameters same 

as in Fig. 4. 
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-0.25 -0.2 -0.15 -0.1 -0.05 0 0.05 0.1 

real part of eigenvalue 

Figure 7: Effect of alternate mistuning on eigenvalues in transonic flow; standard configuration 
No. 10, M = 0.8, /' = 13°, kb = 0.1, N = 12, pitching about mid-chord; structural parameters same 
as in Fig. 4. 
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Figure 8: Comparison of torsional response for a tuned cascade; flat plate geometry, M = 0.7, 
i = 0°, kb = 0.5, N = 12, pitching about leading edge; structural damping ratio = 0.002; structural 
parameters same as in Fig. 4. Lines are calculations from linear theory. 
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Figure 9: Effect of steady loading on torsional response, M =0.7, / = 0°, kb = 0.5, N = 12, 
pitching about leading edge; structural damping ratio = 0.002; structural parameters same as in 
Fig. 4. 



real part of eigenvalue 

Figure A 1: Effect of alternate mistuning on eigenvalues in incompressible flow; flat plate, 
incompressible flow, i = 0°, kb = 0.642, N = 56, pitching about mid-chord; structural parameters 

same as in Fig. 4 (identical with Fig.5 of Ref. 7). 


N AS A/TM— 1 999-209633 


24 






Figure A. 2. Effect of blade alternate mistuning on coupled bending-torsion response in 
incompressible flow; flat plate, incompressible flow, / = 0°, kb = 1.2, N = 56, R = 11, pitching 

about mid-chord; structural damping ratio = 0.002; structural parameters same as in Fig. 4 
(identical with Fig. 1 la of Ref. 7). ~ 



nondimensional forcing frequency 

Figure A. 3: Effect ot blade alternate mistuning on coupled bending-torsion response in 
incompressible flow; flat plate, incompressible flow, / = 0°, kb - 1.2, N = 56, R = 39, pitching 

about mid-chord; structural damping ratio = 0.002; structural parameters same as in Fig. 4 
(identical with Fig. 1 lb of Ref. 7) 
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